I recognize, and fully understand, that this data maybe emotionally difficult to work. My intention is to make these lab relevant, allowing you to gather your own insights directly from new visualizations of the data. Please let me know if you would rather not work with the data.

Learning Objectives

Today’s libaries

library(tidyverse)
library(maps)
library(mapdata)
library(lubridate)
library(viridis)
library(wesanderson)
library(plotly)

Tree Maps

Tree Maps are not geographical maps, but do provide a different perspective on the data. Essentially it is like a pie chart in rectangular form

report_03_11_2020 <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/03-11-2020.csv")) %>%
  rename(Country_Region = "Country/Region", Province_State = "Province/State") %>% 
  group_by(Country_Region) %>% 
  summarise(Confirmed = sum(Confirmed)) %>% 
  mutate(parents = "Confirmed") %>%
  arrange(-Confirmed) 
  
plot_ly(report_03_11_2020,
        type= "treemap",
        values = ~Confirmed,
        labels= ~ Country_Region,
        parents=  ~parents,
        domain = list(column=0),
        name = "Confirmed Cases per Country",
        textinfo="label+value+percent parent")

Creating Country, State and County maps

The ideas for this course tutorial came from multiple examples contributed by Prof. Chris Sutherland and these tutorials

  1. Maps in R using maps by Eric Anderson
  2. geom_maps
  3. Drawing beautiful maps programmatically with R, sf and ggplot2

For our labs two types of approaches will be used to add data to maps. The first that we worked with last week is based on using the longitude and latitude information to add a point at a specific position on a map. The second is two add the information to shapes in a map based on the name of the shape (e.g. states). Although ggmaps is a wonderful tool for mapping using Google Maps and other resources, it is beyond what is needed for now.

As in the previous labs the sources of data are from Github repo for Novel Coronavirus (COVID-19) Cases that supports the dashboard.

For this lab it is important to note that the time series data does not currently have entries for US States. The daily reports include US State and more recently US country/administrative units. Is possible to concatenate the daily reports to create a time series for US States, but cognizant of changes in the formats of the daily reports.

Building maps

Here is a graph containing all the coordinate information. Note this is not summarized by country. Since there are now main points for US counties, there are many points in the US

daily_report <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/04-02-2020.csv")) %>% 
  rename(Long = "Long_") 
    
ggplot(daily_report, aes(x = Long, y = Lat, size = Confirmed/1000)) +
    borders("world", colour = NA, fill = "grey90") +
    theme_bw() +
    geom_point(shape = 21, color='purple', fill='purple', alpha = 0.5) +
    labs(title = 'World COVID-19 Confirmed cases',x = '', y = '',
        size="Cases (x1000))") +
    theme(legend.position = "right") +
    coord_fixed(ratio=1.5)

Zoom in on US 48 states. To do this Alaska, Hawaii and US Territories are filtered. Some US State entries have a Lat and Long of zero, so these are filtered as well.

daily_report <-   read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/04-05-2020.csv")) %>% 
  rename(Long = "Long_") %>% 
  filter(Country_Region == "US") %>% 
  filter (!Province_State %in% c("Alaska","Hawaii", "American Samoa",
                  "Puerto Rico","Northern Mariana Islands", 
                  "Virgin Islands", "Recovered", "Guam", "Grand Princess",
                  "District of Columbia", "Diamond Princess")) %>% 
  filter(Lat > 0)
ggplot(daily_report, aes(x = Long, y = Lat, size = Confirmed)) +
    borders("state", colour = "black", fill = "grey90") +
    theme_bw() +
    geom_point(shape = 21, color='purple', fill='purple', alpha = 0.5) +
    labs(title = 'COVID-19 Confirmed Cases in the US', x = '', y = '',
        size="Cases") +
    theme(legend.position = "right") +
    coord_fixed(ratio=1.5)

Here is a prettier version based on an example by Anisa Dhana

mybreaks <- c(1, 100, 1000, 10000, 10000)
ggplot(daily_report, aes(x = Long, y = Lat, size = Confirmed)) +
    borders("state", colour = "white", fill = "grey90") +
    geom_point(aes(x=Long, y=Lat, size=Confirmed, color=Confirmed),stroke=F, alpha=0.7) +
    scale_size_continuous(name="Cases", trans="log", range=c(1,7), 
                        breaks=mybreaks, labels = c("1-99",
                        "100-999", "1,000-9,999", "10,000-99,999", "50,000+")) +
    scale_color_viridis_c(option="viridis",name="Cases",
                        trans="log", breaks=mybreaks, labels = c("1-99",
                        "100-999", "1,000-9,999", "10,000-99,999", "50,000+"))  +
# Cleaning up the graph
  
  theme_void() + 
    guides( colour = guide_legend()) +
    labs(title = "Anisa Dhana's lagout for COVID-19 Confirmed Cases in the US'") +
    theme(
      legend.position = "bottom",
      text = element_text(color = "#22211d"),
      plot.background = element_rect(fill = "#ffffff", color = NA), 
      panel.background = element_rect(fill = "#ffffff", color = NA), 
      legend.background = element_rect(fill = "#ffffff", color = NA)
    ) +
    coord_fixed(ratio=1.5)

  • Note that in both examples the ggplot funtion borders is used to define the areas in the map

Mapping data to shapes

daily_report <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/04-02-2020.csv")) %>% 
  rename(Long = "Long_") %>% 
  filter(Country_Region == "US") %>% 
  group_by(Province_State) %>% 
  summarize(Confirmed = sum(Confirmed)) %>% 
  mutate(Province_State = tolower(Province_State))
# load the US map data
us <- map_data("state")
# We need to join the us map data with our daily report to make one data frame/tibble
state_join <- left_join(us, daily_report, by = c("region" = "Province_State"))
# plot state map

Using R color palattes

This is a bit of a digression back to Labs 3 and 4, but there are many R color palattes to choose from or you can create your own. In the above a simple gradient is used. The example from Anisa Dhana uses the viridis palatte which is designed to be perceived by viewers with common forms of colour blindness. Here is an example using a different color package - Wes Anderson. …and more Top R Color Palettes to Know for Great Data Visualization

# plot state map
ggplot(data = us, mapping = aes(x = long, y = lat, group = group)) + 
  coord_fixed(1.3) + 
# Add data layer
  geom_polygon(data = state_join, aes(fill = Confirmed), color = "black") +
  scale_fill_gradientn(colours = 
                         wes_palette("Zissou1", 100, type = "continuous"),
                         trans = "log10") +
  labs(title = "COVID-19 Confirmed Cases in the US'")

A look now by the counties using the RColorBrewer

library(RColorBrewer)
# To display only colorblind-friendly brewer palettes, specify the option colorblindFriendly = TRUE as follow:
# display.brewer.all(colorblindFriendly = TRUE)
# Get and format the covid report data
report_03_27_2020 <-   read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/04-02-2020.csv")) %>% 
  rename(Long = "Long_") %>% 
  unite(Key, Admin2, Province_State, sep = ".") %>% 
  group_by(Key) %>% 
  summarize(Confirmed = sum(Confirmed)) %>% 
  mutate(Key = tolower(Key))
# dim(report_03_27_2020)
# get and format the map data
us <- map_data("state")
counties <- map_data("county") %>% 
  unite(Key, subregion, region, sep = ".", remove = FALSE)
# Join the 2 tibbles
state_join <- left_join(counties, report_03_27_2020, by = c("Key"))
# sum(is.na(state_join$Confirmed))
ggplot(data = us, mapping = aes(x = long, y = lat, group = group)) + 
  coord_fixed(1.3) + 
  # Add data layer
  borders("state", colour = "black") +
  geom_polygon(data = state_join, aes(fill = Confirmed)) +
  scale_fill_gradientn(colors = brewer.pal(n = 5, name = "PuRd"),
                       breaks = c(1, 10, 100, 1000, 10000, 100000),
                       trans = "log10", na.value = "White") +
  ggtitle("Number of Confirmed Cases by US County") +
  theme_bw()

If we look at just Massachusetts

daily_report <-   read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/04-02-2020.csv")) %>% 
  rename(Long = "Long_") %>% 
  filter(Province_State == "Massachusetts") %>% 
  group_by(Admin2) %>% 
  summarize(Confirmed = sum(Confirmed)) %>% 
  mutate(Admin2 = tolower(Admin2))
us <- map_data("state")
ma_us <- subset(us, region == "massachusetts")
counties <- map_data("county")
ma_county <- subset(counties, region == "massachusetts")
state_join <- left_join(ma_county, daily_report, by = c("subregion" = "Admin2")) 
# plot state map
ggplot(data = ma_county, mapping = aes(x = long, y = lat, group = group)) + 
  coord_fixed(1.3) + 
# Add data layer
  geom_polygon(data = state_join, aes(fill = Confirmed), color = "white") +
    scale_fill_gradientn(colors = brewer.pal(n = 5, name = "BuGn"),
                         trans = "log10") +
  labs(title = "COVID-19 Confirmed Cases in Massachusetts'")

  • Note the cases on Nantucket and Dukes counties were reported as one value and not included on the graph. There is also an asssigned category that includes 303 Confirmed cases as of 3/31/2020.
daily_report
## # A tibble: 14 x 2
##    Admin2              Confirmed
##    <chr>                   <dbl>
##  1 barnstable                283
##  2 berkshire                 213
##  3 bristol                   424
##  4 dukes and nantucket        12
##  5 essex                    1039
##  6 franklin                   85
##  7 hampden                   546
##  8 hampshire                 102
##  9 middlesex                1870
## 10 norfolk                   938
## 11 plymouth                  621
## 12 suffolk                  1896
## 13 unassigned                270
## 14 worcester                 667

Interactive graphs

In Lab 5 plotly was introduced. It is a great simple way to make interactive graphs with the maps

library(plotly)
ggplotly(
  ggplot(data = ma_county, mapping = aes(x = long, y = lat, group = group)) + 
  coord_fixed(1.3) + 
# Add data layer
  geom_polygon(data = state_join, aes(fill = Confirmed), color = "black") +
    scale_fill_gradientn(colours = 
                         wes_palette("Zissou1", 100, type = "continuous")) +
  ggtitle("COVID-19 Cases in MA") +
# Cleaning up the graph
  labs(x=NULL, y=NULL) +
  theme(panel.border = element_blank()) +
  theme(panel.background = element_blank()) +
  theme(axis.ticks = element_blank()) +
  theme(axis.text = element_blank())
)

Here is an example with the world map

# Read in the daily report
daily_report <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_daily_reports/09-26-2020.csv")) %>% 
  rename(Long = "Long_") %>% 
  group_by(Country_Region) %>% 
  summarize(Confirmed = sum(Confirmed), Deaths = sum(Deaths))
# Read in the world map data
world <- as_tibble(map_data("world"))
# Check to see if there are differences in the naming of countries
setdiff(world$region, daily_report$Country_Region) 
##  [1] "Aruba"                               "Anguilla"                           
##  [3] "American Samoa"                      "Antarctica"                         
##  [5] "French Southern and Antarctic Lands" "Antigua"                            
##  [7] "Barbuda"                             "Saint Barthelemy"                   
##  [9] "Bermuda"                             "Ivory Coast"                        
## [11] "Democratic Republic of the Congo"    "Republic of Congo"                  
## [13] "Cook Islands"                        "Cape Verde"                         
## [15] "Curacao"                             "Cayman Islands"                     
## [17] "Czech Republic"                      "Canary Islands"                     
## [19] "Falkland Islands"                    "Reunion"                            
## [21] "Mayotte"                             "French Guiana"                      
## [23] "Martinique"                          "Guadeloupe"                         
## [25] "Faroe Islands"                       "Micronesia"                         
## [27] "UK"                                  "Guernsey"                           
## [29] "Greenland"                           "Guam"                               
## [31] "Heard Island"                        "Isle of Man"                        
## [33] "Cocos Islands"                       "Christmas Island"                   
## [35] "Chagos Archipelago"                  "Jersey"                             
## [37] "Siachen Glacier"                     "Kiribati"                           
## [39] "Nevis"                               "Saint Kitts"                        
## [41] "South Korea"                         "Saint Martin"                       
## [43] "Marshall Islands"                    "Macedonia"                          
## [45] "Myanmar"                             "Northern Mariana Islands"           
## [47] "Montserrat"                          "New Caledonia"                      
## [49] "Norfolk Island"                      "Niue"                               
## [51] "Bonaire"                             "Sint Eustatius"                     
## [53] "Saba"                                "Nauru"                              
## [55] "Pitcairn Islands"                    "Palau"                              
## [57] "Puerto Rico"                         "North Korea"                        
## [59] "Madeira Islands"                     "Azores"                             
## [61] "Palestine"                           "French Polynesia"                   
## [63] "Western Sahara"                      "South Sandwich Islands"             
## [65] "South Georgia"                       "Saint Helena"                       
## [67] "Ascension Island"                    "Solomon Islands"                    
## [69] "Saint Pierre and Miquelon"           "Swaziland"                          
## [71] "Sint Maarten"                        "Turks and Caicos Islands"           
## [73] "Turkmenistan"                        "Tonga"                              
## [75] "Trinidad"                            "Tobago"                             
## [77] "Taiwan"                              "USA"                                
## [79] "Vatican"                             "Grenadines"                         
## [81] "Saint Vincent"                       "Virgin Islands"                     
## [83] "Vanuatu"                             "Wallis and Futuna"
# Many of these countries are considered states or territories in the JHU covid reports,
# but let's fix a few of them
world <- as_tibble(map_data("world")) %>% 
 mutate(region = str_replace_all(region, c("USA" = "US", "Czech Republic" = "Czechia",  
        "Ivory Coast" = "Cote d'Ivoire", "Democratic Republic of the Congo" = "Congo (Kinshasa)", 
        "Republic of Congo" = "Congo (Brazzaville)")))
# Join the covid report with the map data
country_join <- left_join(world, daily_report, by = c("region" = "Country_Region"))
# Create the graph
ggplotly(
ggplot(data = world, mapping = aes(x = long, y = lat, text = region, group = group)) + 
  coord_fixed(1.3) + 
# Add data layer
  geom_polygon(data = country_join, aes(fill = Deaths), color = "black") +
  scale_fill_gradientn(colours = 
                         wes_palette("Zissou1", 100, type = "continuous")) +
  labs(title = "COVID-19 Deaths")
)

Exercises

  1. For the above graph “COVID-19 Deaths” summarize the counts for each Country on the graph and update the graph to 9/26/2020. You are doing some real life data wrangling. Data is not always in the form that you expected, so it is important to check what the results of each step are. You can summarize the counts for each country and find the median Lat and Long as a way of summarize the Lat and Long from each state. However, the US and several other countries do not have counts. This is because for some US (and other countries) the Lat and Long are NA. One strategies is to simply remove this data (which is fine for this class).
  filter(Lat != "NA")
  filter(Long != "NA")

Using the mean or median(Lat) and (Long) is still not perfect. Some countries are still centered in the ocean. This is ok for ex1. You can use ggplotly to help trouble shoot by putting the Country_Region as text in the hover box

ggplotly(
ggplot(daily_report, aes(x = Long, y = Lat, text = Country_Region, size = Confirmed))...
)

  1. Update Anisa Dhana’s graph layout of the US to 9/26/2020. You may need to adjust the range or change to a linear scale (delete trans=“log”)

  2. Update the above graph “Number of Confirmed Cases by US County” to 9/26/2020 and use a different color scheme or theme

  3. Make an interactive plot using a state of your chosing using a theme different from used in the above exammples.

  4. Create a report with static maps and interactive graphs that is meant to be read by others (e.g. your friends and family). Hide warnings, messages and even the code you used so that it is readable. Included references. Link to the Lab 6 report from your Github site. Submit the link to Moodle. Hide warnings, messages and even the code you used so that it is readable. Included references. Link to the Lab 6 report from your Github site. Submit the link to Moodle. I will talk more about the format on Wed.